Astronomy & Astrophysics manuscript no. Jupiterl l_aa_vlO 


©ESO 2012 


January 24, 2012 





A new vision on giant planet interiors: 
the impact of double diffusive convection 

Jeremy Leconte 1 and Gilles Chabrier 1 ' 

1 Ecole normale superieure de Lyon, CRAL (CNRS), 46 allee d'ltalie, 69007 Lyon, 
Universite de Lyon, France (jeremy.leconte, chabrier @ens-lyon.fr) 

2 School of Physics, University of Exeter, Exeter 

Received 20 February 2009 

Abstract 

While conventional interior models for Jupiter and Saturn are based on the simplistic assumption of a solid core surrounded by 
a homogeneous gaseous envelope, we derive new models with an inhomogeneous distribution of heavy elements, i.e. a gradient 
of composition, within these planets. Such a compositional stratification hampers large scale convection which turns into double- 
diffusive convection, yielding an inner thermal profile which departs from the traditionally assumed adiabatic interior, affecting these 
planet heat content and cooling history. 

To address this problem, we develop an analytical approach of layered double-diffusive convection and apply this formalism to Solar 
System gaseous giant planet interiors. These models satisfy all observational constraints and yield a metal enrichment for our gaseous 
giants up to 30 to 60% larger than previously thought. The models also constrain the size of the convective layers within the planets. 
As the heavy elements tend to be redistributed within the gaseous envelope, the models predict smaller than usual central cores inside 
Saturn and Jupiter, with possibly no core for this latter. 

These models open a new window and raise new challenges on our understanding of the internal structure of giant (solar and extrasolar) 
planets, in particular on the determination of their heavy material content, a key diagnostic for planet formation theories. 
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1. Introduction 

More than 500 planets have now been discovered orbiting stars 
outside our Solar System, spanning a range from a few Earth 
masses to several Jupiter masses. Planets thus seem to be ubiq- 
uitous in nature. These discoveries raise fundamental questions 
about the inner composition, evolution and origin of these bod- 
ies, and about the fundamental properties of planets in general, 
including the ones of our own Solar System. Characterizing their 
internal structure and composition, and from there, better under- 
standing planet formation is one of the major challenges of mod- 
ern astronomy. The determination of the heavy element content, 
for instance, provides key constraints to planet formation mod- 
els, in particular on the efficiency of solid planetesimal accretion 
in the protoplanetary disk to build a planet embryo. While only 
the mass of most of these planets can be derived from observa- 
tions, the mean density can be inferred for about 20% of these 
objects, as they transit their parent star, constraining the planet's 
gross composition (see Baraffe et al. 2010 for a recent review). 
Although providing an important diagnostic, however, this in- 
formation is too limited to determine the element distribution 
within the planet and thus its precise compositional and ther- 
mal structure. Assuming planet formation is a universal process, 
one thus must turn to our Solar System planets, in particular the 
two gas giants, Jupiter and Saturn, which encompass 92% of the 
planetary mass of the Solar System, to derive more detailed in- 
terior models. Indeed, for our own giants, the gravitational mo- 
ments have been determined with high accuracy by the various 
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flyby missions entering their atmosphere and provide stringent 
constraints on their inner element distribution. 

An important question, for instance, is to determine whether 
the heavy elements present in giant planet interiors are located in 
a central core or are mixed in a large fraction into the hydrogen- 
helium (H/He) fluid envelope. In this latter case, a major issue 
is to determine whether convective mixing is sufficiently effi- 
cient to yield a homogeneously mixed envelope or, alternatively, 
if giant planet interiors can exhibit a continuous compositional 
gradient. Exploring such a possibility is crucial to determine (i) 
the maximum amount of heavy elements compatible with obser- 
vational constraints, (ii) the efficiency of heat transport in giant 
planets. These two issues directly impact the planet mechanical 
(density), chemical (composition) and thermal (luminosity, tem- 
perature) structures at a given age, with major consequences on 
our understanding of planet formation and evolution. 

Traditionally, giant planet models have always been based 
on two major specific assumptions concerning their inner pro- 
file, essentially for reasons of simplicity in the planet's model- 
ing (see e.g. Stevenson 1985). It is conventional to assume (i) 
that the inner structure of our giants consists of a few - gener- 
ally 2 to 3 - superposed, well separated, homogeneous regions, 
namely, going from the planet's center to the surface, a central 
solid rocky/icy core, and a surrounding largely dominantly H/He 
gaseous envelope, often split into an inner metallic region and 
an outer atomic/molecular one; the more dense components are 
always supposed to have been accreted first or to have quickly 
settled into the centre under the action of gravity (for planetes- 
imals accreted after the runaway gas accretion); (ii) that plane- 
tary interiors are adiabatic, based on the fact that the giant planet 
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heat flow must be transported by convection (Hubbard 1968). 
All the present determinations of the internal - chemical, me- 
chanical and thermal - structures of the Solar system planets, in- 
cluding their heavy material content, are derived assuming such 
homogeneously stratified, adiabatic interiors (Saumon & Guillot 
2004; Fortney & Nettelmann 2010). 

Giant planet interiors, however, might depart from this con- 
ventional, simplified description, because of complex processes 
for which we lack an accurate description but which may very 
well be at play in real situations (Stevenson 1985). In this pa- 
per, we derive interior models for Jupiter and Saturn which relax 
the aforementioned preconceptions. Instead of the homogeneous 
layer assumption, we explore the possibility of a mixed, inhomo- 
geneous solid-gas interior composition, leading to a heavy mate- 
rial gradient throughout the planet. This in turn tends to suppress 
large scale convection which, due to the double diffusive insta- 
bility (see Sect. 2.2), can turn into either turbulent enhanced dif- 
fusion or layered convection. As both these heat transport mech- 
anisms are fairly inefficient compared to usual convection, this 
compositional gradient thus leads to significant departure from 
global adiabaticity in the interior. As shown below, these mod- 
els do fulfill the planet observational constraints while leading to 
(i) a significantly larger metal content and (ii) significantly larger 
internal temperatures than the one inferred from homogeneously 
stratified adiabatic models. 

This opens a new vision on planet structure, evolution and 
formation efficiency. Such inhomogeneous interior profiles for 
Solar System giant planets had briefly been suggested several 
decades ago by Stevenson (Stevenson 1985) but no attempt 
has ever been made to derive consistent models and to verify 
whether such models would be consistent with the planet var- 
ious observational constraints. This scenario has been revived 
recently in the context of extrasolar planets and has been shown 
to provide a possible or at least complementary explanation for 
the anomalously large observed radii of many of these bodies 
(Chabrier & Baraffe 2007). Indeed, as mentioned above, not 
only an inhomogeneously stratified interior yields a different in- 
terior structure and global metal content, but it decreases heat 
transport efficiency throughout the planet's interior and thus af- 
fects its cooling, thus its mass-radius relationship at a given age, 
a crucial diagnostic to understand (transiting) extrasolar planet 
structure and evolution. Since, as mentioned earlier, only loose 
constraints on the object's internal composition are accessible 
for gaseous exoplanets, it is crucial to verify whether such un- 
conventional internal structures are a viable possibility for our 
own giants. Furthermore, determining the maximum possible 
amount of heavy elements in Jupiter and Saturn and their distri- 
bution within the planet are important diagnostics to understand 
how our own Solar System giants formed. 

In order to address this issue, we first briefly review our cur- 
rent understanding of the double diffusive instability, and of the 
various regimes under which it can occur in Sect. 2. Then, in 
Sect. 3, we derive an analytical formalism, based on a standard 
parametrization of convection, similar to the mixing length for- 
malism, which describes the global transport properties of an 
inhomogeneous convective/diffusive medium. In Sect. 4, we dis- 
cuss the possible existence of an equilibrium size for the convec- 
tive/diffusive layers present in a semi-convective planet and de- 
rive analytical constraints for this equilibrium value. Finally, in 
Sect. 5, we derive semi-convective models of Jupiter and Saturn 
that are consistent with the gravitational moments and with the 
surface abundances measured by the Galileo and Cassini mis- 
sions. This enables us to further constrain the number of possi- 
ble convective/diffusive cells within these planets, and to derive 



new values of the heavy element content and of the core masses 
for our gas giants. 

2. Physical mechanisms leading to inhomogeneous 
density stratification 

2. 1 . Compositional gradient 

In the present calculations, in contrast to all previous planetary 
models, we consider a mass fraction of heavy material Z(m) at a 
depth r{m) within the planet (i.e. at the depth of the iso-density 
surface enclosing a mass m of the planet) that is continuously de- 
creasing from the core to the surface, producing a compositional 
gradient within the gaseous envelope of the planet, 

dlnZ dr 1 dZ H P 

V z = = = -VZ, (1) 

dlnP dlnPZdr Z 

where Hp = -jfes is the pressure scale height. 

The possible origin of such a compositional gradient is an 
open question (Stevenson 1985; Chabrier & Baraffe 2007). In 
the conventional scenario, all the accreted planetesimals are as- 
sumed to directly sink to the core and not to evaporate in the 
envelope, for sake of simplicity. In reality, however, incomplete 
mixing of large planetesimals or dissolution of a substantial frac- 
tion of volatiles and rocks from small solid bodies could occur 
in the envelope during the phase of planetesimal accretion on 
the nascent planet; a substantial amount of ice could as well re- 
main in the envelope (Iaroslavitz & Podolak 2007; Hori & Ikoma 
2011). 

The gradient might also stem from an only partial redistri- 
bution by small scale convective motions of stably layered (sol- 
uble) constituents released by core erosion in the gas-rich en- 
velope during the planet's evolution (Stevenson 1982; Guillot 
et al. 2004), as seems to be supported for water by recent nu- 
merical simulations (Wilson & Militzer 2010). This could be 
enhanced by the immiscibility (phase separation) of an abundant 
enough material (e.g. helium, water) in the dominantly metallic- 
hydrogen envelope (Stevenson & Salpeter 1977a). This would 
change the dynamical properties of (double diffusive) convec- 
tion near the regions of immiscibility and add complexity to the 
problem. Therefore, to avoid extra complication, we only con- 
sider the occurence of this process for soluble material in the 
present paper. 

At last, rapid rotation and/or strong magnetic fields, neces- 
sarily present in Jupiter and Saturn interiors, are known to ham- 
per large-scale convection (Chabrier et al. 2007), possibly lead- 
ing to imperfect mixing of heavy elements in part of the enve- 
lope. 

2.2. Double diffusive convection 

Convective systems in which (rapidly diffusive) heat and (slowly 
diffusive) composition have opposite destabilizing and/or stabi- 
lizing effects tend to develop inhomogeneous density profiles. 
The (de)stabilizing effect of heat can be quantified by the thermal 
gradient, Vj = fg|;, or more precisely by the super adiabaticity 
(Vy - V a( j, where the derivative in V a( j is taken at constant spe- 
cific entropy), and the one of the heavy elements by the gradient 
of mean molecular weight (jj), = J^p. The actual dynamical 
state of the medium depends on the value of the density ratio 

e or Vr-V* 
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Figure 1. Stability diagram for a 
medium presenting a destabilizing 
temperature gradient and a stabilizing 

^_ solute gradient, as a function of in- 

a T Vj — V a( j creasing RZ 1 (decreasing R p ). The usual 
instability Ledoux criterion corresponds 
toR- 1 < 1. 



where a T = - ^\ p _, a u = f^||^ 1 (Stern 1960). 

In the case both the mean molecular weight (due e.g. to a 
higher concentration of salt in salty water) and the temperature 
increase with height (V^ and (V7- - V a d) < 0), the compositional 
gradient is destabilizing while the temperature gradient is stabi- 
lizing. This case is referred to as the fingering case. In that case, 
convective instability develops when R p < 1, which is equiva- 
lent to the Ledoux instability crite rion. But even if 1 < R p < 1/t, 
where r = D/k t is the ratio of solute (D) to thermal (k t ) diffusiv- 
ities, the slower diffusivity of elements compared to heat yields 
the so-called double-diffusive instability which leads to the for- 
mation of salt fingers, and sometimes of thermo-compositional 
staircases, as observed in some parts of the oceans and in labo- 
ratory experiments (Turner 1967). 

The opposite case, referred to as the diffusive case, corre- 
sponds to a fluid exhibiting a positive molecular weight gradient 
(V^ > 0). In that case, the fluid will be convectively unstable if 
this gradient is insufficient to stabilize the system against con- 
vective instability, i.e. if the Ledoux instability criterion 

V r -V ad >^V, o R p > 1, i.e. R p l < 1, (3) 

cej H 

is met 2 . What happens, however, in the regions which are stable 
according to the Ledoux criterion, i.e. R p l > 1, but unstable 
according to the Schwarzschild criterion, 

V r > V ad , (4) 

is less clear, especially in the astrophysical context, where the 
very low values of the Prandlt number (Pr = v/kj, where v is 
the kinematic viscosity) render difficult direct numerical hydro- 
dynamical simulations. In particular, the exact nature of double 
diffusive convection if it occurs - homogeneous oscillatory con- 
vection or layered convection (i.e. uniformly mixed convective 
layers separated by thin diffusive interfaces characterized by a 
steep jump in the mean molecular weight) - remains uncertain. 
Analytical arguments (Radko 2003) and recent 3D hydrodynam- 
ical simulations (Rosenblum et al. 2011; Mirouh et al. 2011), 
however, seem to suggest the picture presented in Fig. 1. When 
the mean molecular weight gradient (oc R p l ) decreases in a sta- 
ble medium, homogeneous oscillatory convection, also called 
turbulent diffusion, first appears for 

Pr + 1 

^min ^ Rp 1 ^ p (oscillatory convection), (5) 

while well defined thermo/compositional layers start to develop 
when 

1 < R~ l < /? m | n (layered convection). (6) 

1 For a perfect gas, = 1 and a T = 1, as will be used hereafter. 

2 As discussed by Rosenblum et al. (201 1), the analogy between the 
fingering and the diffusive case is more apparent when the inverse den- 
sity ratio, R^ 1 , is used, as will be done thoughout the rest of the paper. 



Rj: m corresponds to the point where the solute to heat buoy- 
ancy flux ratio (= y ) stops decreasing when R p 1 increases (see 
Radko 2003 for details). Its exact value, however, depends on 
the characteristics of the medium in a non trivial way and is dif- 
ficult to estimate (Rosenblum et al. 2011; Mirouh et al. 2011). 
For smaller V^, gradients, the medium is unstable according to 
the Ledoux criterion, and the thermal forcing is strong enough 
to force large scale overturning convection. 

Various arguments seem to support, or at least not to exclude, 
the existence of layered convection under planetary conditions 
(Chabrier & Baraffe 2007). Conducting 3D hydrodynamics cal- 
culations over a wide domain of parameter space (Prandlt num- 
ber and atomic to thermal diffusivity ratio) including the regime 
relevant for planetary interiors, Mirouh et al. (201 1) always find 
a domain where y decreases with R~ l , a necessary and suffi- 
cient condition for the layering instability, thus layer formation 
to occur. A central question is then the size of the layers; this 
will be examined in Sect. 4. In any event, both homogeneous 
double-diffusive convection or layered convection - generically 
denominated as "semi-convection" in the following - are found 
to yield thermal and compositional fluxes that are significantly 
smaller than that expected from standard convection. Indeed, the 
presence of diffusive interfaces strongly decreases the efficiency 
of heat transport compared with large-scale, adiabatic convec- 
tion, leading in planet interiors to a significant departure from 
the usual adiabatic profile, as quantified below. 

3. An analytical theory for layered convection 

In order to investigate the impact of such strongly hampered con- 
vection on giant planet internal structure, we developed a simple 
sub-grid model of layered convection. As illustrated on Fig. 2 
and found in simulations (Rosenblum et al. 2011; Mirouh et al. 
2011), we consider that a semi convective zone consists of a 
large number, N\, of well mixed convectively unstable layers of 
size I, separated by thin diffusive interfaces of thermal thick- 
ness 6t, within which the large stabilizing compositional gradi- 
ent completely inhibits convective motions. 

3. 1 . Convective layers 

Within each convective layer, the fluid is expected to follow the 
dynamics found in turbulent Rayleigh-Benard convection within 
a cell of typical height equal to the size of the layer, /. By analogy 
with the mixing length formalism (see details in Appendix A), 
we define a dimensionless mixing length parameter by dividing 
/ by the pressure scale height, a = l/Hp. 

In a laboratory or a numerical experiment, the efficiency of 
the convection is characterized by the thermal Nusselt number, 

A^f^iL, (7) 
F d -Ff 
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Figure 2. Schematic representation of the interiors of Jupiter and Saturn, according to the present study, and of layered convection, 
with the resulting compositional and thermal radial profiles. The abundance of metals is constant within the well mixed convective 
cells of size I, and undergoes a steep variation within the diffusive interfaces of thermal size dp (grey regions). Thanks to this steep 
gradient, these interfaces are stable against convection and energy and matter are transported therein by diffusive processes. Because 
the size of these layers is very small compared with the size of the planet, the mean thermal and compositional gradients ((Vj-) and 
(V„)) can be used in good approximation to infer the global planet structure. 



where by definition the total intrinsic flux (F tot ), the flux trans- 
ported by diffusive processes (Fa), and the diffusive flux that 
would be present in a completely adiabatic zone (F Jj d ) are given 
by (Cox & Giuli 1968; Hansen & Kawaler 1994) 



Fa 



d / 



= Kj- 



Hp 



Vad7 



(8) 



where cp is the heat capacity at constant pressure. It is found 
that, for large Rayleigh numbers, the Nusselt number follows a 
law of the type 



Nu T = C L Ra a *, 



(9) 



where Ra* is a modified Rayleigh number 3 , which is the ratio 
of the strength of the thermal forcing to the one of the radiative 
losses 



Ra* = 



a T gH p 



^(Vr-V, 



a) = (^ 4 /4). 



(10) 



Here, g is the local gravity acceleration, and Np the Brunt- 
Vaisala frequency. 

As convection at very high Rayleigh numbers is difficult to 
study either experimentally or computationally, it is difficult to 
give precise values for a and Cp. For the bounded Rayleigh- 
Benard problem, theoretical models suggest that the exponent of 
the convective flux law, a, could be equal to 1/3 (Garaud et al. 
2010) 4 . However, experiments done by Krishnamurti (1995) 
tend to show that this exponent could be smaller, and as low 



3 Indeed, in the astrophysical context, it is convenient to use Ra* = 
Pr x Ra, where Ra is the usual Rayleigh number. 

4 Although the simulations presented by Rosenblum et al. (2011) 
seem to support exponent values smaller than 1/3, suggesting that in- 
terfaces act as impermeable boundaries, it should be noticed that the 
height of the layers present in their simulations is small compared to a 
pressure scale height. Their conclusion may thus not be valid for larger 
layers. 



as a — 0.2. On the other hand, for homogeneous Rayleigh- 
Benard convection (without boundaries), Garaud et al. (2010) 
showed that the regime predicted by the mixing length theory, 

1/2 

i.e. Nup — Ra* {Cp = 1, a = 1/2; see Appendix A), is recov- 
ered. In the following, we will thus consider 0.2 < a < 0.5 and 
Cl = 1. 

Eq. (9) is sufficient to calculate the flux transported by con- 
vection once the super adiabaticity is known. To compute this 
latter, however, we must first define a quantity which can be 
computed a priori from the local thermodynamical properties of 
the medium and the total internal energy flux to be transported. 
Following Hansen & Kawaler (1994), this convective forcing can 
be defined by 



<5 = Nup xRa*. 



(11) 



Introducing s& = - V a d, we rewrite Eq. (1 1) as <D = On a 4 e d 
where 



On = 



(apgH?^ 



(12) 



It is clear from Eq. (11) that Oo is a local constant of the medium, 
which characterizes its ability to transport energy by convection, 
independently of the mixing length or of the flux to be trans- 
ported (oc g d ). 

Then, from Eqs. (9) and (11), one sees that in a region where 
convection remains efficient enough, 



O = Nup x Ra* = C L Ra* => Ra* 
which yields the super adiabaticity, 

Ep =V T - V ad = 



gd 
Nu T 



C L %a Aa 



!/(!+«) 



l/(l+a) 



(13) 



(14) 



The range of super adiabaticity in the convective layers implied 
by this equation for the various possible exponents a is shown in 



4 



Leconte & Chabrier: A new vision on giant planet interiors 



Fig. 3 (pale red area). As seen, the uncertainty on a leads to a 
large dispersion on this super adiabaticity. In this high convec- 
tive efficiency regime, we can further compute the mean convec- 
tive flux which, by definition, is given by 



pcpT 

Fconv = Kt —r — (Vy 
tip 



(15) 



As mentioned above, this scaling law only applies to the 
vigorous convection regime, i.e O or Rci* » 1. Thus, in 
terms of the layer height, convection remains efficient as long 
asff» a cr i t , where <? cr j t is defined such that Ra+ia^) = 1, and, 
from Eq. (10), the critical layer size is 



( 2 \ 1/4 

= ( k t/n 2 t ) = dPr 



1/4 



(16) 



where d is the lengthscale of the fastest growing mode of the 
linear instability, d = (vk t /N t ) 1/4 (e.g. Baines & Gill 1969). 
Analytical and numerical arguments show that the size of the 
fastest growing layers is equal to 10-100 d > / cl -; t (Radko 2003; 
Rosenblum et al. 2011). The efficient convection regime is thus 
appropriate in the planetary domain, where Pr ~ 10~ 2 - 10~' 
(Chabrier & Baraffe 2007). 

For conditions prevailing in the interior of the actual Jupiter 
(Saturn), the mean thermal diffusivity is kj ~ 5 x 10~ 5 m 2 .s _1 
(Potekhin 1999), O is equal to <D ~ 3 x 10 33 (9xl0 32 ) and 
s d ~ 10 (10; see also Fig. 5), so that a crit = 2 x 10~ 9 (3xl0~ 9 ). In 
the following, all the order of magnitude estimates done through- 
out the text will use these values. 



3.2. Interfaces 

In the interfaces of thermal size St, overturning convection is 
inhibited by the strong jump in molecular weight. However, 
these regions do not need to be in the fully diffusive regime 
(R~ l > (Pr + l)/(Pr + t)), but can also be in the oscillatory 
convection regime characterized by Eq. (5). Indeed, as shown 
by Radko (2005) the condition Rj l > fl m j n , where Rj l is the in- 
verse density ratio within the interface, is a sufficient criterion 
not only to ensure the stability of the interface itself, but also of 
the whole stack of layers which would otherwise merge into one 
homogeneous layer (see Sect. 4.1 for details). 

Therefore, the interface is very likely either in a stable dif- 
fusion state, or in a state of weakly turbulence enhanced diffu- 
sion. Rosenblum et al. (201 1) show that in this regime Nuj < 2, 
meaning that most of the energy is transported by pure diffusion. 
We will thus assume that the thermal gradient to be used in these 
regions is given by the gradient needed to transport the whole 
outgoing energy flux by diffusion: 



1 Hp 
Vd = —F tot . 

Kj pCpl 



(17) 



The thermal diffusivity kj encompasses the contribution of all 
diffusive processes. If diffusion is ensured by photons, as gener- 
ally in most astrophysical objects, the diffusive thermal gradient 
reduces to the so-called radiative gradient, V ra( j (e.g. Hansen 
& Kawaler 1994). However, in the deep interior of giant plan- 
ets and in degenerate bodies, density can be high enough for 
the electrons to become degenerate enough to efficiently conduct 
thermal energy (see Stevenson & Salpeter (1977b) and Chabrier 
& Baraffe (2007) for the characteristic radiative and conductive 
opacities under jovian planet conditions). Along this paper, we 



will use the generic denomination Vd for the diffusive tempera- 
ture gradient, keeping in mind that diffusion is now due to elec- 
tronic or atomic motions, with a characteristic thermal diffusiv- 
ity kj. In the present calculations, we use the conductive thermal 
diffusivities calculated by Potekhin (1999). 

3.3. Mean thermal gradient 

Once we have calculated the thermal gradient in the convective 
zones of size /, and in the diffusive interfaces of size St, we need 
to determine the mean properties of a whole stack of convective- 
diffusive cells. As discussed in Chabrier & Baraffe (2007), since 
the convective plumes must be fed by the diffusive interfaces, 
the thermal convective {N T l ) and diffusive (S t /kt) time scales 
should be similar in each respective layer. Therefore 



M=(«r/PN T r=X"-: /4 =(^c L T^ 



(18) 



Not surprisingly, the size of the interface, St, is related to the 
lengthscale of the most unstable mode of the linear instability by 
St - d/Pr 1 / 4 . In addition, comparing the mean kinetic energy 
of an upwelling eddy with the potential energy barrier created by 
the negative buoyancy in the diffusive interface, we see that the 
above condition also entails that convective overshooting can be 
neglected (see Chabrier & Baraffe 2007 for details). 

Knowing St /I enables us to compute the mean thermal gra- 
dient to be used to compute the planet's structure (see Sect. 5.1) 



St 
l + S 7 



V d + 



/ 



l + S T 



(19) 



Substituting V T and St /I by their expressions in Eqs. (14) and 
(18), and developing the mean gradient to first non vanishing 
order, we find that 



<Vr) * V ad + (V d - V ad ) 



1 



1 



Ra 



; ' 4 + Nu T ) 



+ 



1 



1 



Ra]! 2 ' Nu T RaJ A , 
(20) 



or equivalently 
<v r > - V ad * (V d - V ad ) 



(•/ c J-«- + (<dc£'')- 



(21) 



As expected, two terms appear in the mean super adiabaticity, 
which is shown in Fig. 3 (dark grey region). The first one is due 
to the temperature jumps at each interface and is controlled by 
the convective overturning timescale which determines the size 
of these interfaces. The second one is due to the super adiabatic- 
ity in the convective layers themselves, which is enhanced by 
their small size. The relative contributions of these two terms 
only depend on the properties of the convection, characterized 
by the exponent a. 

Interestingly, if the impact of a on the super adiabaticity in 
the convective layers (pale red area in Fig. 3) is significant, the 
impact on the mean super adiabaticity of the whole stack of lay- 
ers, (St) — ^ad (dark gray region), is rather small. This can be 
understood as follows: for the smallest value of the exponent 
considered here, a = 0.2, convection is very inefficient and most 
of the mean super adiabaticity is contained in the convective lay- 
ers (upper contour of the grey region). In that case, Eq. (19) is 
dominated by the second term, (Vy) ~ Vj-. 

When the exponent is increased, the super adiabaticity in the 
convective layers (oc 0^ a /( 1+a )) drops dramatically, but this en- 
tails a growth of the convective timescale, and thus of the thick- 
ness of the interfaces (see Eq. (18)). For efficient convection, 



5 



Leconte & Chabrier: A new vision on giant planet interiors 





in 




n i 








0.001 


1 

> 


io- 5 




io- 7 




IO" 9 




10 



£d=V d -V ad 
<V r )-V ad 




Figure 3. Range of mean super adiabaticity ((VY) - V ad ) of a 
semi-convective medium as a function of the mixing length pa- 
rameter a = 1 1 Hp, for 0.2 < a < 0.5 (gray area between the 
red dashed curves; O = 3 x 10 33 ; Sd = 10; C L = 1). For 
comparison, the super adiabaticity in a convective (pale red area; 
same parameters) and in a diffusive cell (dotted curve) are also 
shown. As expected a smooth transition between the convective 
and diffusive regimes occurs near a cl -i t (see Appendix A; here 
2 x 10 -9 ^ 



a. 



a ^ 1/4, this effect eventually dominates, and the mean super 
adiabaticity increases again with a. This yields a rather nar- 
row region for the mean super-adiabaticity within the considered 
range of exponent variation. 

Interestingly, to first order, the a = 1/2 and a = 0.2 cases, 
i.e. the two extreme values considered here, yield the same mean 
superadiabaticity. Therefore, the model presented here for semi- 
convection depends only weakly on the precise prescription cho- 
sen for convection. In the following, we will show results ob- 
tained with the mixing length theory, i.e. Cl — 1 and a = 1/2, 
while keeping in mind that for a » a cr j t these results would be 
the same as the one obtained in the a = 0.2 case. 



3.4. Mean solute gradient and flux 

Disregarding convective overshoot, the solute flux (Fz) is deter- 
mined by the transport at each interface. As for heat, the tur- 
bulent transport of solute is expected to be very inefficient, and 
we will assume purely diffusive processes (see Sect. 3.2). The 
impact of a more efficient transport will be discussed in the next 
section. In this regime, the solute flux is thus given by 



-pDVZ/ 



(22) 



where the subscript / describes quantities determined within the 
interface, indicating that the solute flux Fz is determined by the 
element fraction variation rate at each interface. Neglecting the 
small inhomogeneities in the convective layer, the interfacial and 
mean gradients are related by 



l + S z 



<V„> or VZj 



l + 6 z 



<VZ>, 



(23) 



where <VZ) describes the mean gradient of heavy element mass 
fraction, i.e. the value of dZ/dr averaged over several convec- 
tive/diffusive cells, and 5z is the length over which composi- 
tional jump occurs. To constrain the value of 8z, two limit- 
ing arguments can be used. On one side, equating the diffusive 



timescales for heat and solute leads to 

5j V K T 



(24) 



On the other hand, considering that the compositional jump must 
cover the whole interface to ensure its mechanical stability, one 
would expect that 



6t. 



(25) 



The present paper focuses on the impact of double-difffusive 
convection on the planet mechanical and thermal structure. The 
impact on the evolution will be addressed in a forthcoming paper 
(see Sect. 6). Therefore, this uncertainty on the heavy element 
flux, discussed in the next section, does not have a significant 
impact on the results presented in the present study. Determining 
more precisely the solute transport properties in the regime of 
layered convection, however, will be of central importance to 
evolutionary calculations, and 3D hydrodynamical simulations 
in a realistic parameter range are strongly needed (Mirouh et al. 
2011). 

The term which is the relevant quantity when comput- 
ing the equilibrium condition of the interface, is given by 



l + S ZjJ dlnp 
S z oZ 



<vz>. 



(26) 



P.T 



The precise value of 



d ln/j 

~3z~ 



P.T 



depends on the precise chemical 



composition of the considered heavy element, but is typically 
around unity; this is the value we will use in numerical estimates 
below. We will also use D m 10~ 7 - 10~ 8 m 2 .s _1 (i.e. t w 10~ 2 - 
10" 1 ), appropriate for giant planet interior conditions Stevenson 
& Salpeter 1977b). 

4. Layer size: analytical point of view 

4.1. Existence of an equilibrium height 

The determination of the size of the layers, when layered con- 
vection is occurring, is a complex task. The problem is rendered 
even more difficult by the fact that the small layers that initially 
form tend to merge into larger layers. The question is then to 
know if and when layer merging ends. 

In the fingering case, i.e. when the solute is the destabiliz- 
ing actor, Radko (2005) derived a criterion for the linear stability 
of thermo-compositional staircases against merging. In this pic- 
ture, an initially inhomogeneous medium starts forming layers 
if its density ratio (R p ) is smaller than the density ratio, R m \ n , 
for which the heat to solute buoyancy flux ratio (= y) stops de- 
creasing when R p increases. Then, layers are unstable and merge 
as long as the density ratio within the interface between them 
(Rj) is smaller than R m j„. Then, as Rj increases with the layer 
height (under some conditions), it eventually reaches R m i a and 
the merging process stops. At this stage, an equilibrium height 
is reached by the convective layers and the staircase is fully equi- 
librated. 

The linear stability analysis of Radko (2005) can be redone 
in the diffusive case, by simply accounting for the fact that the 
signs of the various density and temperature differences and 
of the fluxes must be changed. Then, after some lengthy but 
straightforward algebra, one can show that layers are unstable 
and merge when 



dy 



< o RJ 1 < R~ 



(27) 
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which defines /£? . As for the fingering case, this linear stability 
analysis shows that, if layers form, their merging will stop when 
their height reaches a finite equilibrium value. 

4.2. Global constraints 

Although precisely estimating this equilibrium value is difficult, 
some strong limits on the layer size can be derived theoretically. 
Let us consider a stack of layers extending over a zone of size 
L, and define a global gradient (VZ) * -AZ/L, where AZ is 
the difference between the mass ratios of heavy element at the 
bottom and at the top of the semi-convective zone. 

On one hand, the mean molecular weight gradient in all the 
interfaces, V^j, must be high enough to satisfy the stability cri- 
terion discussed in sections 2.2 and 3.2. This implies 



aj 



V>tf ra ! n x(v d -v ad ). 



(28) 



Substituting V^j by its expression in Eq. (26) and 6z by using 
Eqs. (24) and (18), the criterion (28) for a planetary scale (L w 
R p ) semi-convective zone reads 



mm 



aj 



d /rUVd-Vad) 1 



*r (d> /C L )^AZ ^ 



(29) 



P.T 



where a m i n denotes a lower limit for the layer size. Under the 
present conditions in the interiors of our gas giants (see end of 
Sect. 3.1 and Table 1), this yields a min « 2 x 1(T 9 x (/T] n ) 1+a , 
with the less restrictive constraint being obtained for R 1 =1, 

° mm ' 

as summarized in Table 1 . However, as the existence of layers 
allows /T! n to be as high as (1 + Pr)/(r + Pr) m 10 1 - 10 2 , these 
constraints could be severely tightened, as showed in Table 1. 
This simply confirms that, in order for layered convection to 
be stable, convective cells must remain larger than the diffusive 
interfaces, and the medium is always in the convective regime, 
a » a crit . 

Table 1. Numerical constraints on the layer height in the two 
limiting cases of convection (a = 0.2 and 0.5; Cx = 1) for the 
following conditions (representative of Jupiter interior): <Dq 



3 x 10 33 , e d = 10, k t 



5x 10- 5 m 2 .s-\ 



D 



5 x lO-W.s -1 . 





constraint 






a = 0.2 


a = 0.5 


^min 


stability 


R-J. 

nun 

R-J 

nun 


= 1 

_ l+Pr 
r+Pr 


2.4x10^ 
6.0xl0~ 7 


2.4x10-* 
2.4xl0- 6 


^max 


homogeneization 




= 1 

= 3 


1.0x10 4 
2.5xl0~ 5 


1.4xl0- 3 
2.6xl0" 4 


Q'min 


observational 
constraints* 






Jupiter 
Saturn 


3 x 10" 5 

4 x 10~ 6 



* See Sect. 5 

On the other hand, the solute gradient within the planet will 
be homogenized within a typical timescale 



tz 



pAZR f 
\Fz\ 



Using Eq.(22), (23) and (24), and taking <VZ) 
Eq. (30) becomes 



R 1 , 



t z ~ 



Rldz 

D I \C L 



-l/O+a) 



(30) 
AZ/fl p , 

(31) 



Therefore, in order to avoid a complete homogenization of giant 
planet interiors in less than 5 Gyr, their present age, a must be 
smaller than a max » 1 x 10~ 4 for a = 0.2 and 1.4 x 10~ 3 for 
a = 0.5. Note however that, as mentioned in Rosenblum et al. 
(2011), turbulent transport due to the double-diffusive instabil- 
ity can yield compositional Nusselt numbers, Nu^, around 2- 
4, yielding even stronger constraints, as summarized in Table 1. 
These small values of a justify a posteriori the approximation 
of continuous thermal and heavy element profiles when consid- 
ering the planet's entire internal structure. Note, however, that 
layered inhomogeneities could be dynamically regenerated over 
time. In that case, layered convection will be an ongoing process 
in the planet's interior. 

Considering the possibility that the ratio of the compositional 
to the thermal size of the interface does not scale as the square 
root of the ratio of the diffusivities, and therefore using Eq. (25) 
instead of Eq. (24), yields an increase of both a m m and a max by 
a factor yfcrjD 1 ^' * 10. 

This analysis shows that, for the age of the Solar System and 
for the conditions prevailing in gas giant interiors, there exists 
a range of layer sizes for which ongoing layered convection is 
a viable mechanism. According to our estimate, this range is 
relatively large and corresponds to : 



10 y - 10~ 6 < a < 10~ 4 - 10" 



(32) 



The uncertainty on the lower and upper bounds are respectively 
due to our poor knowledge of the behavior of R m \ n , and of the 
solute flux at low Prandlt number. Even given these uncertain- 
ties, the fact that convection is always in a regime of relatively 
high Rayleigh number appears to be robust prediction. To make 
an attempt to overcome these limitations, we examine in Sect. 5 
how observational data can narrow this possible domain of a, 
by further constraining the degree of super adiabaticity, and thus 
the size of the layers, in our Solar System gas giants. 

4.3. An alternative scenario 

We stress that the aforementioned constraints apply only to lay- 
ered convection, and does not preclude the possibility that, un- 
der some conditions, double diffusive convection may manifest 
itself under the form of homogeneous double diffusive convec- 
tion and act like a turbulent diffusion (Rosenblum et al. 2011). 
However, for this to happen, the criterion (28) must be verified 
in the inefficient convection regime, St » I (see Appendix A). 
The following criterion must then hold 

-H P (VZ) * AZ ^ > R-J m x (V d - V ad ). (33) 

As AZ < 1 and R 1 > 1, by definition, a zone of turbu- 

mm ' J ' 

lent diffusion cannot extend over the entire planet's scale unless 
V d < V ad , i.e. if the whole object is nearly diffusive in the first 
place. 

In an object with a heat flux high enough to be convectively 
unstable, the semi-convective zone must then be confined to a 
fraction of the planet, in particular, but not necessarily, near an 
immiscibility region or a phase transition for instance. In this 
case, the total size of the zone must verify L/R p < e d ' (~ 1/10 
in Jupiter), condition for which a large enough jump in the heavy 
element mass fraction can be sufficient to stabilize the whole 
zone against convection and open a diffusive buffer in the inte- 
rior (where * V d ). From the global point of view of the 



5 Analogously to the thermal Nusselt number, the flux of heavy ele- 
ments, for a given Nu^, is given by F z = - Nu^ x pDVZ. 
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planet, this would act as a composition, temperature and entropy 
nearly discontinuity. This possibility for the existence of such a 
diffusive-like buffer in the interior of our gas giants, as a conse- 
quence of double-diffusive instability, must be kept in mind. 

In the following, however, we will not consider this scenario 
any further, and we will only consider the effect of a planetary 
scale layered-convection zone. 

5. Numerical results for Solar System giant planets 

In this section, we examine whether the presence of semi- 
convection in Jupiter and Saturn interiors can be consistent with 
the various available observational constraints. We first derive 
homogeneous reference interior models in Sect. 5.2. Then, in 
Sect. 5.3, we incorporate our model for layered convection into 
the standard method used to compute interior structure models of 
rotating gaseous planets (presented in Sect. 5.1) and determine 
the area of the composition/layer size space parameter which 
is consistent with observed gravitational moments and surface 
abundances. 

5.1. Hydrostatic equilibrium and figures of the planet 

Solar System giant planets are rapidly rotating bodies (the pe- 
riod of rotation is about 10 hours), with the centrifugal potential 
representing about 10% of the gravitational potential. This mod- 
ifies the hydrostatic equilibrium condition between the pressure 
gradient and the gravitational force in the interior, which now 
writes 



VP = -pV(V c + V rot ), 



where 



and 



V G (r) 



= -G f ^-dV 
J |r~r'| 

Vrot(r,e) = - \ cutoWf 
Jo 



(34) 



(35) 



(36) 



denote respectively the gravitational and centrifugal potentials, 
with differential rotation « p (£), where £ is the distance from the 
position r to the rotation axis, and G the gravitational constant. 
In the present study, o> p is assumed to be constant and given 
by the magnetospheric rotation rate. Because of the symmetry 
of the centrifugal potential with respect to both the rotation axis 
and the equatorial plane, surfaces of equal densities for these 
objects are supposed to be generalized ellipsoids of revolution 
whose exact shape is given by 



r(r, 6) = r 



1 + 



S2nW^2n(cOS0) 



where r is the mean radius of the equipotential, P2„ are the usual 
Legendre polynomials, 9 is the colatitude and the S2n are a set 
of figure functions. These latter can be derived using the theory 
of figures for rotating bodies (Zharkov & Trubitsyn 1978), and 
must be solved iteratively with the set of perturbed ID hydro- 
static equilibrium equations 



dP 


1 Gm 


dm 


~ ~4n^~ 


df 


1 


dm 


4nf 2 p ' 


dT 


t dp„ 




= -^v r 


dm 


P dm 



6tt? 



(38) 



(39) 



(40) 



where m is the mass enclosed in the equipotential of mean radius 
r, w p is the rotation rate of the planet, <p w (r) is a second order 
correction due to the centrifugal potential, which depends on the 
figure functions. As discussed in Sect. 3, the prescription to be 
used for V7- is determined by the energy transport processes. 

The departure from sphericity of the iso-density surfaces re- 
sults in a perturbation of the external gravity field Voir, 8) that 
writes 

V O (r,0) = - — jl-Xh? J2nP2n(cOS0)}, (41) 
n=l * ' 



p(r, 6) r 2n P 2 n(cos ff) d 3 r, (42) 



where r is the radial distance from the center of the planet, M the 
mass of the planet, R eq the equatorial radius, 9 the colatitude, Pi n 
are Legendre polynomials of order 2n and denote the gravita- 
tional moments, that can be computed once the figure equations 
have been solved. The measured gravity moments provide strin- 
gent constraints on the density profile and the possible layering 
within these planets. 

As, in practice, Legendre polynomial expansions are trun- 
cated at a given order n, a closure equation is provided by the 
equation of state (EOS) of the mixture along the planet's interior 
profile. Such an EOS is generally given by the so-called ideal 
volume law for the mixture: 



1 X Y Z 

- = — + 1 , 

P PX PY PZ 



(43) 



where X, Y and Z denote the mass fractions of H, He and heavy 
elements, respectively. For the H/He fluid, the most widely 
used EOS is the Saumon-Chabrier-vanHorn EOS (Saumon et al. 
1995; SCvH). For the heavy material, we have used the "Rock" 
EOS of Hubbard & Marley (1989) for silicates and the "Ice" 
ANEOS equation of state (Thompson & Lauson 1972) for 
volatiles (CH4, NH3, H2O). The impact of the differences be- 
tween various EOS's on exoplanet structure and evolution has 
been explored in Baraffe et al. (2008). 

Once such EOS's, P[p(Xj)], are specified, structure models 
with various compositions are calculated by solving iteratively 
the aforementioned hydrostatic equilibrium condition for a ro- 
tating body and the third-order level-surface theory (Zharkov 
& Trubitsyn 1978) to obtain a model which reproduces the ob- 
served values of the radius, R eq , and gravitational moments Ji 
and J4 measured by the Pioneer and Voyager missions (see Table 
2). 



(37) 5.2. Reference homogeneous models 



In conventional giant planet models, the abundances of heavy 
elements are chosen to be constant in the gaseous H/He enve- 
lope, with a possible discontinuity at the transition between the 
molecular and metallic regions (Chabrier et al. 1992). Under 
the actual conditions found in Jupiter and Saturn, the thermal 
gradient that would be needed to transport the whole flux by dif- 
fusive processes, Vj, is always larger than V a d and the whole 
interior is convective in the homogeneous case according to the 
Schwarzschild criterion. As convection is very efficient (see 
Sect. 3), the super adiabaticity needed to transport the outgoing 
energy is on the order of 10~ 8 - 10~ 9 , so that the structure can be 
solved by setting V7- = V a( j in Eq. (40). 

In order to have a reference case, we use the formalism de- 
scribed in Sect. 5.1 to derive homogeneous, adiabatic interior 
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Table 2. Observed characteristics of Solar System gaseous gi- 
ants (Guillot 2005 and references therein; the numbers in paren- 
theses are the uncertainty in the last digits of the value). 





Jupiter 


Saturn 


M p [KPkg] 


18.986112(15) 


5.684640(30) 


R eq [10 7 m] 


7.1492(4) 


6.0268(4) 


tf P oi[10 7 m] 


6.6854(10) 


5.4364(10) 


Ptat [10 4 s] 


3.57297(41) 


3.83577(47) 


Tlbm [K] 


165.(5) 


135.(5) 


F tot [W.m- 2 ] 


5.44(43) 


2.01(14) 


J 2 x 10 2 


1.4697(1) 


1.6332(10) 


J 4 x 10 4 


-5.84(5) 


-9.19(40) 




2-4 


2-8 


(y/(x + y)) atm 


0.238(50) 


0.215(35) 



models representative of the usual 2-layer composition. As we 
use the interpolated SCvH EOS, we do not consider the effect 
of a Plasma Phase Transition, and we are thus left with only 
two free parameters, namely the core mass (M c ) and the metal 
mass fraction in the gaseous envelope (Z env ) 6 . The temperature, 
density and pressure profiles of our best representative homoge- 
neous models of Jupiter and Saturn are shown in Fig. 4 (solid 
curves). These are composed of a solid core of mass M c = 3.9 
and 25.6 M m surrounded by a H/He gaseous envelope with a con- 
stant metal fraction Z env = 0.11 and Z env = 0.05 for Jupiter and 
Saturn, respectively (these results are summarized in Table 3). 
These reference models yield interior enrichment that are consis- 
tent with previous determinations (Chabrier et al. 1992; Saumon 
& Guillot 2004; Guillot 2005). 

5.3. Inhomogeneous models 

We now derive semi-convective, inhomogeneous interior mod- 
els for Jupiter and Saturn. We stress that all these models are 
consistent, within the observational uncertainties, with the mea- 
sured gravitational moments of Jupiter and Saturn (see Table 2; 
Campbell & Synnott 1985; Campbell & Anderson 1989). 

An additional constraint on the outermost value of the com- 
positional gradient is provided by the surface abundance of 
heavy elements in the planets measured by the 1995 Galileo 
Entry Probe mission. Indeed, elemental abundances of the atmo- 
spheres of solar giant planets are observed to differ significantly 
from each other and from the solar composition, being enriched 
by a factor ~ 2 - 4 and ~ 2 - 8 with respect to the Sun's atmo- 
sphere for Jupiter and Saturn, respectively, as shown in Table 2 
(Guillot 2005). Moreover, the planet's total mean abundances of 
H and He (X and Y) must recover the values of the protosolar 
nebula, i.e. Y/(X + f) * 0.275. 

In our calculations, the adjustable parameters to fulfill all 
these constraints are chosen to be the mass of the core (M c ), 
the mean heavy element mass fraction in the gaseous envelope 
(Zenv), an d the global compositional variation in the envelope 
(AZ env , the difference between the metal mass fraction just above 
the central core and the one in the atmosphere). To assess the 
robustness of our results with respect to the equation of state 
chosen to describe the thermodynamics of the heavy material, 



6 As already found by Chabrier et al. (1992) and Saumon & Guillot 
(2004), with the SCvH EOS, simple homogeneous models such as our 
reference case cannot reproduce 74 to better than a few percents error. 




m/M p 



Figure 4. Pressure (black), temperature (red) and density (blue) 
profiles as a function of depth (expressed by the Lagrangian co- 
ordinate, i.e. the mass m), for the reference adiabatic (solid 
curves) and semi-convective (dashed curves) cases for Saturn (a) 
and Jupiter (b). The increased thermal gradient due to the inef- 
ficient heat transport in the semi-convective case (with N\ = 10 4 
for Jupiter and 10 4 5 for Saturn) strongly increases the internal 
temperature. This causes a partial redistribution of the core ma- 
terial within the gaseous envelope. 

we derived several sets of models for which the composition of 
the core varies from pure ice to pure rock. 

5.3.1. Number of layers 

From a macroscopic point of view, an important quantity 
describing layered convection is the number of convective- 
diffusive layers, N\. This number is roughly equal to the ratio of 
the size of the semi-convective zone, comparable to the planet's 
radius, R p , if this zone extends over the whole planet, to the 
height of a typical convective/diffusive cell, I + 6t- As shown 
in Sect. 4, in the regime of interest, I + 6j ~ I, and N\ ~ R v /l. 
Because Hp w R p in the deep interior, the number of layers in 
the planet is thus approximately equal to N\ ~ or -1 , and in the 
following we will always refer indifferently to either a or 

N\ = a' 1 = Hp/l. (44) 

As shown in Sect. 4, N\ is constrained to lie within the range 
10 2 ~ 4 < N\ < 10 6 ~ 9 . Note that, given the small size of the 
diffusive-convective layers compared with the size of the planet, 
the discontinuous (staircase-like) temperature and composition 
profiles can be well approximated by continuous mean thermal 
and compositional gradients ((Vj-) and (V^), respectively) to de- 
termine the planet's global structure, as illustrated in Fig. 2. 

This possible range of numbers of layers is further con- 
strained by our numerical calculations, which show that, in order 
to reproduce our giant planet observational constraints, no more 
than ~ 2.5 x 10 5 layers (a m ; n ~ 4x 10~ 6 ) can in reality be present 
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Figure 5. Conductive (top red curves) and adiabatic (bottom 
black curves) thermal gradient profiles for Jupiter and Saturn. 
The dashed curves correspond to the mean thermal gradient for 
models with 10 4 5 and 10 5 4 layers for Jupiter and Saturn, respec- 
tively. The shaded area represents the allowed range of super 
adiabaticity in presence of semi-convection, consistent with the 
observational constraints. 



in Saturn and ~ 3 x 10 4 (a m i n 3 x 10~ 5 ) in Jupiter (see Table 1). 
Indeed, a larger number of layers leads to so high temperatures 
in the interior that the induced mean density decrease can not be 
counterbalanced by an increase of the heavy element mass frac- 
tion compatible with the observed surface abundances. This is 
due to the fact that, the larger the number of layers, the smaller 
the size of each convective cell, reducing the maximum height a 
convective eddy can travel to transport heat before being stopped 
by the negative buoyancy present in the diffusive interface. 

A large number of layers thus decreases convective heat (and 
composition) transport efficiency. This leads to an increase of 
the mean super adiabaticity, as portrayed on Fig. 5, which in turn 
immediately implies a rise of the internal temperature, as illus- 
trated on Fig. 4. It is important to stress that super adiabaticity 
is the physical quantity most directly constrained by the data. 
Thus, whereas the allowed range of number of layers (or equiv- 
alently of layer sizes) may slightly depend on the model used 
to parametrize semi-convection (see Sect. 3), the allowed range 
of super adiabaticity displayed in Fig. 5 should remain weakly 
affected. 

5.3.2. Enhanced heavy material enrichment 

The pressure, density and thermal profiles obtained in the most 
extreme semi-convective case compatible with the observational 
constraints discussed above are shown in Fig. 4 (dashed curves). 
As seen on the figure, and as expected from the above discus- 
sion, the non-adiabatic envelope profile obtained in the semi- 
convective case yields substantially higher internal temperatures 
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Figure 6. Metal abundance profiles as a function of depth (ex- 
pressed by the Lagrangian mass m) for Saturn (a) and Jupiter 
(b) for different numbers of layers. The abundance increases 
with the number of layers, to keep the density profile unchanged 
when convection becomes inefficient. The extreme cases (10 5 4 
layers for Saturn and 10 4 5 layers for Jupiter) correspond to the 
semi-convective profiles portrayed in Fig. 4. 



than the usual adiabatic assumption, as heat and material redis- 
tributions are partly inhibited by diffusive processes. The pres- 
sure and density profiles, on the other hand, remain barely af- 
fected, being strongly constrained by the gravitational moments. 
Hence, at basically fixed density profile, a higher temperature 
profile must be compensated by a larger amount of heavy mate- 
rial within the envelope. This is illustrated in Fig. 6, where we 
show the abundance profiles, as calculated in Appendix B, cor- 
responding to semi-convective models with different numbers of 
layers. The bottom curve (solid) in each panel corresponds to 
models with 1000 layers while the other curves correspond to a 
gradually increasing number of layers. 

Therefore, in order to compensate the radius increase (den- 
sity decrease) due to the hotter interior, semi-convection yields 
a significantly larger total metal content compared with conven- 
tional homogeneous models. This can be seen in Fig. 7, which 
shows the amount of heavy elements in the core and envelope 
for the various cases discussed here, as summarized in Table 
3. For Saturn, up to 50 M @ of heavy elements could be present 
in the planet while for Jupiter the heavy material content could 
reach 63 M e . This corresponds to about 25 and 10 times the so- 
lar abundances, respectively 7 . Since these values only depend 
on the allowed amount of super adiabaticity, they should not 
strongly depend on our modeling of diffusive/convective trans- 



7 Note that the abundances of heavy elements brought to Jupiter and 
Saturn, in particular water, could already be enriched compared with 
the solar value (Gautier et al. 2001) 
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port, as mentioned above. In contrast, the maximum amount 
of heavy elements compatible with the observational constraints 
for the homogeneous, adiabatic models, is about 30M ffl for 
Saturn and 40 M ffl for Jupiter, in agreement with previous studies 
(Saumon & Guillot 2004). 

But semi-convection does not only increase the global metal 
content, it also yields a completely different distribution of heavy 
elements. While the global enrichment of the planet is increased 
in the inhomogeneous models, the mass of the central core is 
decreased, as heavy elements are preferentially redistributed in 
the gaseous envelope. In the case of Saturn, the vertical spread 
in core mass at fixed number of layers illustrated in Fig. 7 is ob- 
tained when varying the core composition from pure ice (top) 
to pure rock (bottom). In Jupiter the inferred core mass is too 
small for the equation of state to make a significant difference. 
One could wonder why the homogeneous case is not continu- 
ously recovered when a tends toward 1. This slightly counter 
intuitive effect is due to the fact that, at least when using the 
SCvH EOS, completely homogeneous models (central core plus 
a fully homogeneous envelope) cannot in general reproduce both 
the observed Ji and J4 (Chabrier et al. 1992; Saumon & Guillot 
2004). Thus, if we relax the constant Z condition in the enve- 
lope, the presence of a compositional gradient and of a smaller 
core appears to be the best solution to reproduce observational 
data, even in the absence of any additional super-adiabaticity. 

For Jupiter, models can be found that match the gravitational 
moments without the presence of a central, completely differen- 
tiated core (red dots on the bottom right of Fig. 7). Such cases 
yield an atmospheric metallicity Z atm ~ 4 - 5Z G . The fact that 
the possible erosion of the core mass would have been more effi- 
cient in Jupiter than in Saturn might stem from the larger energy 
flux available in Jupiter (Guillot et al. 2004). 

Table 3. Heavy element content for Jupiter and Saturn inferred 
from the various models consistent with these constraints within 
the quoted observational uncertainties. 
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Figure 7. Mass range of heavy elements in the core (M c ) and 
in the envelope (Mz env ) consistent with all observational con- 
straints, for different numbers of layers, for Jupiter (bottom 
right) and Saturn (upper left). The open dots at the upper left of 
each region correspond to the homogeneous interior models. As 
the number of semi-convective layers increases, the efficiency of 
convection decreases, and the heavy element mass fraction in- 
creases to counteract the radius increase induced by the planet's 
higher internal temperature. The metals initially present in the 
core are then redistributed within the envelope. For Jupiter, so- 
lutions with no core at all (M c = 0) can be found for the non 
adiabatic models (red dots). 



Jupiter 



Saturn 



Region Amount of heavy elements (M e ) 

Homogeneous model 
Envelope 36 4.7 
Core 3.9 25.6 

Total 40 30.3 

Semi-convective models 
Envelope 41-63.5 10-36 
Core 0-0.5 10-21 

Total 41-63 26-50 



6. Prospect for giant planet evolution 

While the aim of present study is to examine and to constrain 
the properties of transport by semi-convection inside Jupiter and 
Saturn at the present time, the impact of non-adiabatic interi- 
ors on the cooling of these planets, and of giant planets in gen- 
eral, remains to be explored. This requires more cumbersome 
evolutionary calculations, for which initial conditions will be of 
crucial importance, as will be explored in a forthcoming study. 
Note that evolution will add an additional constraint, namely that 
the planet cooling rates yield the correct properties at the age of 
the solar system, which will put more stringent constaints on the 
range of possible layer sizes/numbers. Without going into such 



detailed calculations, however, the following points can be men- 
tioned. 



6. 1 . Merging of the layers 

As mentioned in Sect. 4, soon after they form, layers are ex- 
pected to merge, leading to larger layers, until the layer height 
reaches either the planet size, yielding a standard adiabatic in- 
terior, or an equilibrium value, which is the semi-convective 
case considered here. Numerical simulations by Radko (2005) 
and Rosenblum et al. (2011) tend to show that the equilibra- 
tion timescale of the staircase is much shorter than the typical 
timescale for the evolution of the planet. 

Thus, if an equilibrium height is reached, as discussed in 
Sect. 4.1, this rather quick equilibration timescale suggests that 
the layer size should remain roughly constant during the evolu- 
tion, or change slowly with the mean properties of the medium 
(e.g. the mean molecular weight gradient). A precise prescrip- 
tion for the height of the layers being yet lacking and demand- 
ing a more precise knowledge of the heat transport properties of 
layered convection under astrophysical conditions, it seems rea- 
sonable, as a first guess, to use a constant size, whose value has 
been constrained in the afore sections, throughout the evolution. 
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6.2. Initial heavy element distribution 

Conventional models based on fully adiabatic thermal profiles 
notably lead to cooling times about 15% longer than the age of 
the Solar System for Jupiter (Fortney et al. 2011). In princi- 
ple, the hotter non-adiabatic internal structures suggested in the 
present paper will prolong the cooling and thus worsen the prob- 
lem. 

However, in the case of the erosion of an initially large core, 
part of the gravitational work will be spent eroding the core 
and mixing the material upward and will thus not contribute to 
the total luminosity, quickening the cooling. All these effects 
must be properly accounted for to infer the appropriate cooling 
timescale. 

In addition, if Jupiter and Saturn initial cores were allowed 
to be relatively large (> 10 M e ), the corresponding high surface 
density of solids in the protosolar nebula will quicken the forma- 
tion timescale in the conventional core accretion scenario, help- 
ing solving the related formation timescale problem (Pollack 
et al. 1996). Finally, since, in the present scenario, some of 
the ablated material from the accreted planetesimals during the 
planet's early formation stages remains distributed throughout 
the envelope, this will (i) reduce the heating due to gravitational 
energy release produced by the infalling planetesimals on the 
planet embryo and (ii) increase the envelope mean molecular 
weight. Both effects will cause the protoplanet to contract more 
quickly, shortening again the planet's formation timescale in the 
conventional core accretion scenario (Pollack et al. 1996). A cor- 
rect exploration of the impact of inhomogeneous interiors upon 
giant planet history thus necessitates to investigate the conse- 
quences not only on the thermal evolution but also on the forma- 
tion process. 

7. Conclusion and perspective 

In this paper, we have first developed an analytical approach of 
layered convection, based on a standard parametrization of con- 
vection similar to the MLT formalism. This formalism allows 
a quantitative determination of the expected number of diffu- 
sive layers, or equivalently of the average characteristic mixing- 
length parameter, in a semi-convective planet interior character- 
ized by a given total flux and a given thermal (and composi- 
tional) diffusivity. Furthermore, this formalism allows an ex- 
act determination of the characteristic thermal gradient in the 
presence of double-diffusive convection, and thus of the related 
amount of super-adiabaticity within the planet's interior. 

Using this formalism, we have computed semi-convective in- 
terior models of Jupiter and Saturn. We have shown that a strat- 
ified internal structure for Solar System gaseous giants, with a 
compositional gradient of heavy material extending over a sub- 
stantial fraction of the planet, is a viable hypothesis, as such 
models can fulfill all the observed gravitational and atmospheric 
constraints for these planets. This new possibility differs from 
the conventional description of giant planet interiors, assumed 
to be composed of 2 main superposed, well identified layers of 
homogeneously distributed material, namely a solid core sur- 
rounded by a dominantly gaseous H/He envelope. The conse- 
quences of the present giant planet interior description are mul- 
tiple. Namely, 

- (i) our jovian planets might be significantly more enriched in 
heavy elements than previously thought, 

- (ii) their interior temperature, thus heat content, might be 
much larger than usually assumed, 



- (iii) the inner temperature profile could significantly depart 
from the usually assumed adiabatic profile. 

We stress that these conclusions do not depend on the precise 
model used to describe double diffusive convection. Besides 
directly affecting our conventional vision of giant planet me- 
chanical, compositional and thermal structures, these results 
have profound impacts on our understanding of planet formation 
and cooling properties. Indeed, the revised possible maximum 
amount of heavy material bears direct consequences on the de- 
termination of the efficiency of solid planetesimal accretion dur- 
ing planet formation in the protoplanetary nebula, suggesting an 
early and efficient capture of planetesimals for our, and probably 
extrasolar as well, giant planets. Moreover, the larger heat con- 
tent and the departure from adiabaticity, as well as the possibility 
of significant core erosion from an initially large core, directly 
impact the planet cooling histories. Departure from adiabaticity, 
in particular, implies less efficient heat transport, a direct con- 
sequence of the inhibited convective motions due to a persistent 
compositional gradient, and thus a smaller heat flux output rate 
than assumed in the conventional approach. 

These results open a new window, and raise new challenges, 
on our present understanding of planet structure, formation and 
evolution. Importantly, the viability of such stratified interior 
models for our Solar System gas giants directly applies to the 
case of extrasolar planets, reinforcing the possibility that such a 
lower heat flux output could at least partly explain the anoma- 
lously large radius of several transiting "hot Jupiters" (Chabrier 
& Baraffe 2007). Indeed, it seems that invoking an extra source 
of (tidal, kinetic or magnetic) energy dissipation in these object 
interiors can not completely solve this "radius anomaly" puzzle 
and that an alternative or complementary process is necessary 
(Laughlin et al. 2011). Unconventional, inhomogeneous non- 
adiabatic planetary interiors, as suggested in the present study, 
might provide the missing piece of the puzzle. 
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Appendix A: The case of the mixing length theory 

Here, we briefly discuss the particular case of the Mixing Length 
Theory (MLT; Hansen & Kawaler 1994). This case can be re- 
covered in the efficient convection regime by the more general 
model presented in Sect. 3, by setting Cl — 1 and a = 1/2, but 
we show below that the MLT formalism can be extended to the 
inefficient convection regime. 

Hansen & Kawaler (1994) showed that in the fr « 1 
regime, recasting their Eq. (5.60) in the notation of Sect. 3, the 
modified Rayleigh number verifies 

<D = Nu T x Ra* = Ra+ + {Raj 1 ■ o-(R ai ,)f , (A.l) 



with 



o-iRa*) = — — — ( Vl +4/to* - l) , 



(A.2) 



where ctNt is the growth rate of a convective eddy and the in- 
verse of the convective time. In general, these equations can be 
solved numerically. It can be easily verified that in the limit of 
efficient convection, i.e. Ra+ » 1, o~ — > 1, meaning that ra- 
diative losses are negligible and that the convective time tends 
toward Nt- In this limit, a simple expansion of Eq. (A.l) yields 
<J> w Raj 1 and Nu T * Raj 1 as expected in the standard MLT 
formalism (Hansen & Kawaler 1994). 

In the 1 » Ra+ » Pr regime, however, developing Eq. (A.2) 

1 /2 

yields cr * Raj . Then, because the convective time is dramati- 
cally increased, Eq. (18) rewrites in that case 

(Sr/j} * (*r/ z 2 ^ Nt )- V2 * Ra-: !1 » 1. (A.3) 

The size of the diffusive interfaces thus grows until the convec- 
tive layers eventually turn into a completely diffusive medium 
for which (VY) = Vd. 



can be present at the same time and either compete or contribute 
constructively. 

In our model, we consider an ideal mixture of heavy el- 
ements within a H/He envelope whose H/He mass ratio is 
kept constant and equal to its value in the protosolar nebula, 
(H/He) proto . This implies 



X + Y + Z= 1, and 



X+Y \X + Y 



0.275, (B.l) 



proto 



everywhere in the planet's gaseous envelope. We are then left 
with only one degree of freedom. Following previous calcula- 
tions (Chabrier et al. 1992), for sake of simplicity and in order 
to have a flexible determination of the metal enrichment and a 
thermodynamically consistent EOS in the gaseous phase, we ap- 
proximate the metal mass fraction by an effective helium mass 
fraction (1") in the H/He EOS. For the core, the metal mass frac- 
tion is correctly described by the appropriate water and silicate 
EOS mentioned in the text. 

The various element mass fractions, then the corresponding 
metal enrichment, are thus inferred from the relation 



1 



1 -Z 



p(P,T,Y') Pwm (P,T) p z (P,T)' 



(B.2) 



which gives Z at each depth along a given model P-T profile 
(Chabrier et al. 1992). The hydrogen and helium mass fractions 
are then derived using 



Y = (l-Z) 



X+Y 



proto 



and 



X = 1 - Y - Z. 



(B.3) 



(B.4) 



In this simple model, a Z gradient thus necessarily yields a 
competing inhomogeneous helium distribution within the planet. 
Because the mean molecular weight of a H/He mixture at fixed 
temperature and pressure only depends on Y/{X + Y), only the 
Z variations need to be considered to compute in our sim- 
plified model. The Z(m) profile is then integrated to obtain the 
total amount of heavy elements mixed in the gaseous layers for 
each planetary model, as portrayed in Fig. 7. In the most general 
case, with an intrinsic inhomogeneity of the helium distribution, 
caused for instance by its immiscibility in metallic hydrogen, 
both the Y and Z gradients would have to be properly calculated. 



Appendix B: Computation of the heavy element 
content 

The mean molecular weight gradient needed to stabilize the fluid 
against large scale convection can be caused by an inhomoge- 
neous distribution of both helium (Y(m)) and metals (Z(w)) in 
the hydrogen (X(m))-hch medium. In practice, both gradients 
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